% 定义函数 f(x)
f = @(x) 1./(1 + x.^2);

% 导数 f'(x)
df = @(x) -2*x./(1 + x.^2).^2;

% Theorem 3.57 的插值点
t1 = -6:1:4;

% Corollary 3.58 的插值点
t2 = (-5:1:4) - 0.5;

% 计算插值点上的函数值和导数值
y1 = f(t1);
dy1 = df(t1);
y2 = f(t2);

% Theorem 3.57: Complete Cubic Cardinal B-splines
% 注意：这里简化了实现，具体的 B-spline 构建方式可能需要根据定理进行调整
cs1 = csape(t1, [dy1(1) y1 dy1(end)], 'clamped');

% Corollary 3.58: Quadratic Cardinal B-splines
% 注意：这里简化了实现，具体的 B-spline 构建方式可能需要根据推论进行调整
cs2 = csape(t2, [y2(1) y2 y2(end)], 'second');

% 定义误差计算点
x_points = [-3.5, -3, -0.5, 0, 0.5, 3, 3.5];

% 计算原函数在这些点的值
f_values = f(x_points);

% 计算 Cubic B-spline 在这些点的值
cubic_spline_values = fnval(cs1, x_points);

% 计算 Quadratic B-spline 在这些点的值
quadratic_spline_values = fnval(cs2, x_points);

% 计算插值误差 E_S(x)
E_cubic = abs(cubic_spline_values - f_values);
E_quadratic = abs(quadratic_spline_values - f_values);

% 绘制两个 B-spline 误差图在同一个图中
figure;
plot(x_points, E_cubic, '-o', 'LineWidth', 2, 'DisplayName', 'Cubic B-spline Error');
hold on; % 保持当前图，以便在同一图上绘制下一组数据
plot(x_points, E_quadratic, '-x', 'LineWidth', 2, 'DisplayName', 'Quadratic B-spline Error');
hold off;

% 添加图例、标题和轴标签
legend('show');
title('Cubic vs. Quadratic B-spline Interpolation Error');
xlabel('x');
ylabel('Error');
grid on;
